 #!/usr/bin/perl -w
use strict;

use Bio::DB::SeqFeature::Store;
my $promotor_size =500;

my $dir = shift @ARGV || die "Need a directory of genomes";
$dir =~ s/\/$//;
opendir(DIR, $dir) || die $!;
my @genomes;
for my $file ( readdir(DIR) ) {
    next if $file =~ /^\./; # skip dot files
    push @genomes, $file;
}

for my $genome ( @genomes ) {
    warn("processing $genome ($dir/$genome)\n");
    my $dbh = Bio::DB::SeqFeature::Store->new(-adaptor=>'berkeleydb',
					      -dir    => "$dir/$genome",
					      );
    # iterate through all the TEs
    my $iterator = $dbh->get_seq_stream(-type => ['match_part:repeatmasker']);
    while( my $TE_feature = $iterator->next_seq ) {
	
	# check to see what this overlaps in gene space	
	my $segment = $dbh->segment($TE_feature->seq_id, 
				    $TE_feature->start =>  $TE_feature->end);
	my %overlaps;
	if( my @mRNA_features = $segment->features(-type => ['mRNA']) ){
	    for my $mRNA ( @mRNA_features ) {
		$overlaps{'mRNA'}++;

		for my $exon ( $mRNA->get_SeqFeatures() ) {
		    if( $exon ) {
			if( $exon->overlaps($TE_feature) ) {
			    push @{$overlaps{exon}}, $exon;
			}
		    }
		}
	    }
	    if( exists $overlaps{'exon'} )  {
		print  "$TE_feature ->  overlaps exon\n";
	    }
	}
	# check to see if it overlaps introns or exons?

	# check to see if it overlaps a gene promotor by expanding
	# the search location
	$segment = $dbh->segment($TE_feature->seq_id, 
				 $TE_feature->start - $promotor_size,
				 $TE_feature->end   + $promotor_size);
	if( my @genes = $segment->features(-type => ['gene']) ) {
	    for my $gene ( @genes ) {
		push @{$overlaps{'gene_prom'}}, $gene;
	    }
	}
	if( exists $overlaps{'gene_prom'} && ! exists $overlaps{'mRNA'} )  {
	    print "$TE_feature ->  overlaps promotor\n";
	}
    }
}
